#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
泰坦尼克生存预测 - 机器学习经典入门项目
重点：完整的数据分析、特征工程、模型训练流程
"""

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
from sklearn.model_selection import train_test_split, cross_val_score, GridSearchCV
from sklearn.preprocessing import StandardScaler, LabelEncoder
from sklearn.ensemble import RandomForestClassifier, GradientBoostingClassifier
from sklearn.linear_model import LogisticRegression
from sklearn.svm import SVC
from sklearn.metrics import accuracy_score, precision_score, recall_score, f1_score
from sklearn.metrics import confusion_matrix, classification_report, roc_auc_score, roc_curve
import warnings
warnings.filterwarnings('ignore')

plt.rcParams['font.sans-serif'] = ['Microsoft YaHei']
plt.rcParams['axes.unicode_minus'] = False

print("🚢 泰坦尼克生存预测项目开始")
print("=" * 60)

# 1. 数据加载与探索
print("1. 数据加载与探索")
print("-" * 40)

# 从网络加载泰坦尼克数据集（使用seaborn内置数据集）
try:
    df = sns.load_dataset('titanic')
    print("✅ 成功加载泰坦尼克数据集")
except:
    print("⚠️ 无法加载内置数据集，使用模拟数据")
    # 创建模拟数据
    np.random.seed(42)
    n_samples = 891
    df = pd.DataFrame({
        'survived': np.random.choice([0, 1], n_samples, p=[0.62, 0.38]),
        'pclass': np.random.choice([1, 2, 3], n_samples, p=[0.24, 0.21, 0.55]),
        'sex': np.random.choice(['male', 'female'], n_samples, p=[0.65, 0.35]),
        'age': np.random.normal(29, 14, n_samples).clip(0, 80),
        'sibsp': np.random.poisson(0.5, n_samples).clip(0, 8),
        'parch': np.random.poisson(0.4, n_samples).clip(0, 6),
        'fare': np.random.lognormal(2.5, 1.2, n_samples).clip(0, 512),
        'embarked': np.random.choice(['S', 'C', 'Q'], n_samples, p=[0.72, 0.19, 0.09])
    })

print(f"数据集形状: {df.shape}")
print(f"\n数据前5行:")
print(df.head())

print(f"\n数据基本信息:")
print(df.info())

print(f"\n缺失值统计:")
print(df.isnull().sum())

print(f"\n目标变量分布:")
survival_rate = df['survived'].mean()
print(f"生存率: {survival_rate:.2%} ({df['survived'].sum()}/{len(df)})")

# 2. 数据可视化分析
print("\n2. 数据可视化分析")
print("-" * 40)

fig = plt.figure(figsize=(15, 12))

# 子图1: 生存率分布
plt.subplot(2, 3, 1)
sns.countplot(data=df, x='survived', palette=['red', 'green'])
plt.title('生存情况分布')
plt.xlabel('是否生存 (0=否, 1=是)')
plt.ylabel('人数')

# 子图2: 性别与生存关系
plt.subplot(2, 3, 2)
sns.countplot(data=df, x='sex', hue='survived', palette=['red', 'green'])
plt.title('性别与生存关系')
plt.xlabel('性别')
plt.ylabel('人数')

# 子图3: 船舱等级与生存关系
plt.subplot(2, 3, 3)
sns.countplot(data=df, x='pclass', hue='survived', palette=['red', 'green'])
plt.title('船舱等级与生存关系')
plt.xlabel('船舱等级')
plt.ylabel('人数')

# 子图4: 年龄分布
plt.subplot(2, 3, 4)
sns.histplot(data=df, x='age', hue='survived', kde=True, palette=['red', 'green'])
plt.title('年龄分布与生存关系')
plt.xlabel('年龄')
plt.ylabel('密度')

# 子图5: 票价与生存关系
plt.subplot(2, 3, 5)
sns.boxplot(data=df, x='survived', y='fare', palette=['red', 'green'])
plt.title('票价与生存关系')
plt.xlabel('是否生存')
plt.ylabel('票价')

# 子图6: 热力图 - 特征相关性
plt.subplot(2, 3, 6)
correlation_matrix = df.select_dtypes(include=[np.number]).corr()
sns.heatmap(correlation_matrix, annot=True, cmap='coolwarm', center=0)
plt.title('特征相关性热力图')

plt.tight_layout()
plt.savefig('titanic_eda.png', dpi=300, bbox_inches='tight')
plt.show()

# 3. 特征工程
print("\n3. 特征工程")
print("-" * 40)

# 复制原始数据
X = df.copy()

# 处理缺失值
print("处理缺失值...")
if 'age' in X.columns and X['age'].isnull().sum() > 0:
    X['age'].fillna(X['age'].median(), inplace=True)
if 'embarked' in X.columns and X['embarked'].isnull().sum() > 0:
    X['embarked'].fillna(X['embarked'].mode()[0], inplace=True)

# 创建新特征
print("创建新特征...")
# 家庭大小
X['family_size'] = X.get('sibsp', 0) + X.get('parch', 0) + 1
# 是否独自旅行
X['is_alone'] = (X['family_size'] == 1).astype(int)
# 年龄分组
X['age_group'] = pd.cut(X['age'], bins=[0, 12, 18, 35, 60, 100], 
                        labels=['儿童', '青少年', '青年', '中年', '老年'])
# 票价分组
X['fare_group'] = pd.qcut(X['fare'], q=4, labels=['低', '中低', '中高', '高'])

# 编码分类变量
print("编码分类变量...")
label_encoders = {}
for col in ['sex', 'embarked', 'age_group', 'fare_group']:
    if col in X.columns:
        le = LabelEncoder()
        X[col] = le.fit_transform(X[col].astype(str))
        label_encoders[col] = le

# 选择特征
features = ['pclass', 'sex', 'age', 'sibsp', 'parch', 'fare', 'embarked', 
           'family_size', 'is_alone', 'age_group', 'fare_group']
features = [f for f in features if f in X.columns]

X_features = X[features]
y = X['survived']

print(f"最终特征数量: {len(features)}")
print(f"特征列表: {features}")

# 4. 数据预处理
print("\n4. 数据预处理")
print("-" * 40)

# 划分训练集和测试集
X_train, X_test, y_train, y_test = train_test_split(
    X_features, y, test_size=0.2, random_state=42, stratify=y
)

print(f"训练集大小: {X_train.shape}")
print(f"测试集大小: {X_test.shape}")
print(f"训练集生存率: {y_train.mean():.2%}")
print(f"测试集生存率: {y_test.mean():.2%}")

# 数据标准化
scaler = StandardScaler()
X_train_scaled = scaler.fit_transform(X_train)
X_test_scaled = scaler.transform(X_test)

# 5. 模型训练与评估
print("\n5. 模型训练与评估")
print("-" * 40)

# 定义多个模型
models = {
    '逻辑回归': LogisticRegression(random_state=42),
    '随机森林': RandomForestClassifier(n_estimators=100, random_state=42),
    '梯度提升': GradientBoostingClassifier(random_state=42),
    '支持向量机': SVC(probability=True, random_state=42)
}

results = {}

for name, model in models.items():
    print(f"\n{name}模型训练中...")
    
    # 训练模型
    if name in ['逻辑回归', '支持向量机']:
        model.fit(X_train_scaled, y_train)
        y_pred = model.predict(X_test_scaled)
        y_proba = model.predict_proba(X_test_scaled)[:, 1]
    else:
        model.fit(X_train, y_train)
        y_pred = model.predict(X_test)
        y_proba = model.predict_proba(X_test)[:, 1]
    
    # 计算评估指标
    accuracy = accuracy_score(y_test, y_pred)
    precision = precision_score(y_test, y_pred)
    recall = recall_score(y_test, y_pred)
    f1 = f1_score(y_test, y_pred)
    auc = roc_auc_score(y_test, y_proba)
    
    results[name] = {
        'accuracy': accuracy,
        'precision': precision,
        'recall': recall,
        'f1': f1,
        'auc': auc,
        'model': model,
        'y_pred': y_pred,
        'y_proba': y_proba
    }
    
    print(f"准确率: {accuracy:.3f}")
    print(f"精确率: {precision:.3f}")
    print(f"召回率: {recall:.3f}")
    print(f"F1分数: {f1:.3f}")
    print(f"AUC分数: {auc:.3f}")

# 6. 模型比较与可视化
print("\n6. 模型比较与可视化")
print("-" * 40)

fig = plt.figure(figsize=(15, 10))

# 子图1: 模型性能对比
plt.subplot(2, 3, 1)
metrics = ['accuracy', 'precision', 'recall', 'f1', 'auc']
metric_names = ['准确率', '精确率', '召回率', 'F1分数', 'AUC']
colors = ['skyblue', 'lightgreen', 'lightcoral', 'gold', 'lightblue']

x = np.arange(len(models))
width = 0.15

for i, metric in enumerate(metrics):
    values = [results[name][metric] for name in models.keys()]
    plt.bar(x + i*width, values, width, label=metric_names[i], color=colors[i])

plt.xlabel('模型')
plt.ylabel('分数')
plt.title('模型性能对比')
plt.xticks(x + width*2, models.keys(), rotation=45)
plt.legend()
plt.grid(True, alpha=0.3)

# 子图2: 最佳模型ROC曲线
best_model_name = max(results.keys(), key=lambda x: results[x]['auc'])
best_result = results[best_model_name]

fpr, tpr, _ = roc_curve(y_test, best_result['y_proba'])
plt.subplot(2, 3, 2)
plt.plot(fpr, tpr, color='darkorange', lw=2, label=f'ROC曲线 (AUC = {best_result["auc"]:.3f})')
plt.plot([0, 1], [0, 1], color='navy', lw=2, linestyle='--', label='随机猜测')
plt.xlabel('假正率')
plt.ylabel('真正率')
plt.title(f'{best_model_name} ROC曲线')
plt.legend()
plt.grid(True, alpha=0.3)

# 子图3: 最佳模型混淆矩阵
plt.subplot(2, 3, 3)
cm = confusion_matrix(y_test, best_result['y_pred'])
sns.heatmap(cm, annot=True, fmt='d', cmap='Blues', 
           xticklabels=['未生存', '生存'], yticklabels=['未生存', '生存'])
plt.xlabel('预测标签')
plt.ylabel('真实标签')
plt.title(f'{best_model_name}混淆矩阵')

# 子图4: 特征重要性（随机森林）
plt.subplot(2, 3, 4)
if hasattr(models['随机森林'], 'feature_importances_'):
    importances = models['随机森林'].feature_importances_
    indices = np.argsort(importances)[::-1]
    
    plt.barh(range(len(indices)), importances[indices], color='lightblue')
    plt.yticks(range(len(indices)), [features[i] for i in indices])
    plt.xlabel('特征重要性')
    plt.title('随机森林特征重要性')

# 子图5: 预测概率分布
plt.subplot(2, 3, 5)
plt.hist([results['随机森林']['y_proba'][y_test == 0], 
          results['随机森林']['y_proba'][y_test == 1]], 
         bins=20, alpha=0.7, label=['未生存', '生存'], color=['red', 'green'])
plt.xlabel('预测概率')
plt.ylabel('频数')
plt.title('预测概率分布')
plt.legend()

# 子图6: 交叉验证结果
plt.subplot(2, 3, 6)
rf_model = RandomForestClassifier(n_estimators=100, random_state=42)
cv_scores = cross_val_score(rf_model, X_features, y, cv=5, scoring='accuracy')
plt.boxplot(cv_scores)
plt.ylabel('准确率')
plt.title('5折交叉验证结果')
plt.grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig('titanic_model_comparison.png', dpi=300, bbox_inches='tight')
plt.show()

# 7. 超参数调优
print("\n7. 超参数调优")
print("-" * 40)

# 随机森林参数调优
param_grid = {
    'n_estimators': [50, 100, 200],
    'max_depth': [3, 5, 10, None],
    'min_samples_split': [2, 5, 10]
}

rf = RandomForestClassifier(random_state=42)
grid_search = GridSearchCV(rf, param_grid, cv=5, scoring='accuracy', n_jobs=-1)
grid_search.fit(X_train, y_train)

print(f"最佳参数: {grid_search.best_params_}")
print(f"最佳交叉验证分数: {grid_search.best_score_:.3f}")

# 使用最佳参数重新训练
best_rf = grid_search.best_estimator_
y_pred_best = best_rf.predict(X_test)
accuracy_best = accuracy_score(y_test, y_pred_best)
print(f"调优后测试集准确率: {accuracy_best:.3f}")

# 8. 业务洞察与结论
print("\n8. 业务洞察与结论")
print("-" * 40)

print("🔍 关键发现:")
print(f"1. 最佳模型: {best_model_name} (AUC: {best_result['auc']:.3f})")
print(f"2. 测试集准确率: {best_result['accuracy']:.3f}")
print(f"3. 最重要的特征: {features[np.argmax(best_rf.feature_importances_)]}")

print("\n💡 生存预测洞察:")
print("- 女性生存率显著高于男性")
print("- 头等舱乘客生存机会更大")
print("- 儿童和青少年有更高的生存概率")
print("- 家庭规模适中的乘客生存率更高")

# 9. 项目总结
print("\n9. 项目总结")
print("-" * 40)

print("🎯 本项目实现了完整的机器学习流程:")
print("1. 数据探索与可视化")
print("2. 特征工程与数据预处理")
print("3. 多种模型训练与比较")
print("4. 超参数调优与模型优化")
print("5. 业务洞察与结果分析")

print("\n📊 技术要点:")
print("- 处理分类变量和缺失值")
print("- 创建有意义的衍生特征")
print("- 使用多种评估指标全面评估模型")
print("- 交叉验证确保模型稳定性")

print("\n🚀 下一步改进方向:")
print("- 尝试更复杂的特征工程")
print("- 集成学习方法（如投票分类器）")
print("- 使用更高级的模型（如XGBoost）")
print("- 部署为Web应用")

print("\n" + "=" * 60)
print("✅ 泰坦尼克生存预测项目完成!")
print("🎉 恭喜你完成了第一个完整的机器学习项目!")

print("\n📚 推荐下一步项目:")
print("- 鸢尾花分类（多分类问题）")
print("- 房价预测（回归问题）")
print("- 手写数字识别（计算机视觉入门）")

print("\n💡 运行结果说明:")
print("- 查看生成的PNG图片了解数据分析和模型性能")
print("- 重点关注特征重要性和业务洞察")
print("- 思考如何将学到的知识应用到其他项目")